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ABSTRACT 


Betavoltaic nuclear batteries offer a promising alternative energy source that harnesses the power 
of beta particles emitted by radioisotopes. To satisfy the power demands of microelectromechanical 
systems (MEMS), 3D structures have been proposed as a potential solution. Accordingly, this paper 
introduces a novel 3D °Ni-SiC-based P*PNN* structure with a multi-groove design, avoiding the 
need for PN junctions on the inner surface, and thus reducing leakage current and power losses. 
Monte Carlo simulations were performed considering the fully coupled physical model to extend 
the electron-hole pair generation rate to a 3D structure, enabling the efficient design and 
development of betavoltaic batteries with complex 3D structures. As a result, the proposed model 
produces the significantly higher maximum output power density of 19.74 W/cm? and 
corresponding short-circuit current, open-circuit voltage, and conversion efficiency of 8.57 uA/cm?, 
2.45 V, and 4.58%, respectively, compared with conventional planar batteries. From analysis of the 
carrier transport and collection characteristics using the COMSOL Multiphysics code, we provide 
deep insights regarding power increase, and elucidate the discrepancies between the ideal and 
simulated performances of betavoltaic batteries. Our work offers a promising approach for the 
design and optimization of high-output betavoltaic nuclear batteries with a unique 3D design, and 


serves as a valuable reference for future device fabrication. 
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1. Introduction 


With remarkable performance, compact size, and high-volume production, 
microelectromechanical systems (MEMS) have attracted tremendous interest for various 
applications ranging from mobile electronics to deep-sea surveys, implantable medical devices, and 
autonomous wireless sensor networks. However, the development of micro-batteries capable of 
outputting power of 1-100 uW, long service life, and fitting within the size range of 1 um to 10 mm 
has become a critical challenge for MEMS applications [1-3]. Recently, betavoltaic nuclear batteries 
have emerged as a highly attractive energy option for MEMS applications, with the potential 
advantages of long operational life, high energy density, ultra-miniature size, and strong anti- 


interference [4-7]. 


The performance of betavoltaic nuclear batteries is governed by radioisotope source 
characteristics, device geometry, and semiconductor converter properties. The corresponding 


relationship can be mathematically expressed as: 
QVocFF 
Pn = AE avens(1 =r) = AEavgs (1 -r) = (1) 


where A is radioisotope activity determined by the radioisotope half-life, Eavg is the average B-decay 
energy dictating the total input power driving the converter, ys is the fraction of B-energy emitted to 
the total decay energy of the radioisotope source, r is the reflectivity coefficient, y: is the energy 
conversion efficiency of the semiconductor converter, Q is the carrier collection efficiency, Voc is 
the open-circuit voltage, FF is the fill factor, and e is the effective ionization energy to generate an 
electron-hole pair (EHP). The overall conversion efficiency can be expressed as: mtomns(1-r)yc. The 
output power of betavoltaic nuclear batteries depends on significant losses, such as self-absorption 
loss (1-75) in the radioisotope source, backscattering loss (r) at the semiconductor interface, and 
carrier collection loss (1-Q) outside the depletion region [8-10]. To maximize the energy coupling 
between the radioisotope and semiconductor, it is crucial to increase their face-to-face area while 
reducing self-absorption and backscattering effects. Additionally, 7, depends on the material 
properties of the semiconductor, including bandgap, depletion region width, and carrier diffusion 
length. Q, Voc, and FF can be improved via optimization of doping concentrations and junction 
width to increase the depletion region width and carrier diffusion length, ensuring that the B-particle 


penetration depth matches the converter scale length. 


As the conventional planar configuration of betavoltaic nuclear batteries uses only one side 
of the radioisotope source attached to the semiconductor converter, the output power and conversion 
efficiency are limited. A high activity (4) or input power (Pin) requires a thick source layer, but this 
results in stronger self-absorption effects and smaller ys, saturating A-7s. The directional loss is 
approximately 50%, and the backscattered loss significantly reduces the device efficiency by up to 


25% [10, 11]. Although efforts have been made to improve the conversion efficiency and output 


power of betavoltaic nuclear batteries, including use of a reflector to reduce directional loss and 
backscattered loss [12], optimization of the junction depth and doping concentration to increase the 
carrier collection efficiency [13-16], use of an extra graded N layer to reduce radiation-induced 
EHPs recombination loss [17], and adoption of radioisotope sources with higher particle energy to 
increase input power [18]; their reported values remain limited by the effective loading activity of 
the radioisotope source as well as the coupling efficiency between the source and device, i.e., the 
limitation of A7,(1-r). As a result, the respective output power is only 0.1—50 nW for practical and 
tested batteries [18-22] and 10-400 nW for theoretically predicted ones [13-15, 19, 23, 24], falling 
short of meeting the power requirements of MEMS [25]. 


To achieve higher performance in betavoltaic nuclear batteries, 47s(1-r) can be increased over 
a wide range using specific types of radioisotope sources and semiconductor materials. Increasing 
the specific surface area of the converter can enable a higher loading amount of radioisotope sources, 
leading to larger Ay.(1-r) and higher output power. Compared with two-dimensional (2D) planar 
structures, three-dimensional (3D) structures with a larger specific surface area can significantly 
increase the output power density owing to three factors: (i) more radioisotope sources can fill the 
interspace of 3D structures, (ii) the thinning of radioisotope sources in 3D structures significantly 
reduces the self-absorption effects, and (iii) the combination of a radioisotope source with 3D 
structures leads to the interaction of beta particles with the converter in all directions, increasing the 


collection efficiency of the beta particles. 


In recent years, the use of 3D structures in betavoltaic nuclear batteries has demonstrated 
significant potential for improving the corresponding specific surface area and conversion efficiency; 
thus, they are promising options for meeting MEMS power demands [26-29]. However, 
conventional 3D structures require preparation of PN or PIN junctions on the inner surface of the 
microstructure, which significantly affects the leakage current and output performance [28]. This is 
a predominant reason why the device performance is far from ideal, even up to several orders of 
magnitude, so full use of the 3D structure is challenging. Moreover, while Monte Carlo simulations 
[15, 30] and empirical formulas [31, 32] are widely employed to calculate the distribution of the 
EHP generation rate [G(x)] in betavoltaic nuclear batteries with 2D diode structures, the EHP 
generation rate is rarely evaluated in 3D structured converters combined with radioisotope sources 
distributed in 3D space. Most currently available models for G(x) in 3D structured converters only 
describe specific structures with fixed source and device geometries, and cannot accurately evaluate 
the EHP distribution in 3D structures, highlighting the need for a precise model to advance the 
development of 3D batteries. 

This paper introduces a novel approach that addresses these issues, utilizing a °Ni-SiC-based 
(P*PNN*) structure with a multi-groove design, enabling the epitaxial growth of graded P and N 


layers on the substrate without the need of preparing PN junctions on the inner surface of the 


microstructure. This approach has the potential to significantly reduce leakage current and power 
losses, thereby narrowing the gap between theoretical predictions and experimental results. In 
addition, a novel formulaic model is proposed for calculating the complex EHP generation rate in 
3D-structured betavoltaic nuclear batteries. The model considers the intricate 3D structure of the 
converter and radioisotope source, enabling accurate evaluation of the EHP distribution in all 
possible 3D structures resulting from changes in their geometries. Our fully coupled model, 
combined with the COMSOL Multiphysics code, involves the entire physical process of carrier 
evolution, including B-particle generation, energy deposition, radiation-carrier generation, drift- 
diffusion, and recombination. This study provides valuable insight into the internal mechanisms of 
carrier transport, collection characteristics, and power increase. From the results, our approach 
demonstrates maximized output power density with optimized source thickness, converter geometry, 
doping concentration, and width of each region, outperforming conventional planar batteries. Our 
novel G(x) model provides a critical tool for designing and optimizing 3D structured betavoltaic 


nuclear batteries, with potential applications in other betavoltaic nuclear batteries. 
2. Model and method 
2.1 Device structure 
Fig. 1(a) depicts a °°Ni-SiC-based betavoltaic nuclear battery with a multi-groove structure, 
characterized by the ridge width [d], ridge spacing [¢], and groove depth [H_source]. The converter 


comprises four layers: a P*- SiC layer, graded P-SiC layer, graded N-SiC layer, and N*- SiC layer, 
with thicknesses of H_P*, H_P, H_N, H_N' (the N*-layer thickness is slightly larger than that of 


H_N’), respectively. The “Ni source was filled in multi-grooves, surrounded by converters at the 
sides and bottom (front and rear sides not shown), and enclosed by a metal electrode at the top. This 
design reduces the directional and electrode shielding losses compared with conventional planar 
diode structures. The rectangular top section of the device has the area of 1 cm x 1 cm. To prevent 
PN junction shorting and metal-semiconductor contact formation, techniques such as nitride 
passivation, plasma-enhanced chemical vapor deposition, and atomic layer deposition can be 
employed to grow a thin insulating layer (e.g., Siz3Na, SiO2, and Al203) with a thickness of 10-100 
nm [33-36]. This layer effectively blocks the flow of electrons or holes, achieving insulation 
between the source and semiconductor devices, while having a negligible impact on the energy 


deposition of the source decay energy in the device. 
To ensure good Ohmic contact, the doping concentrations of heavily doped P*- SiC and N*- 


SiC layers are 10'° and 10!8 cm”, respectively. The P- and N-SiC layers serve as the core regions of 


the betavoltaic nuclear battery, generating an internal electric field [2] to separate the radiation- 
I 


induced EHPs. These layers are lightly doped to obtain larger depletion region width [Wa] and 
minority diffusion length [Zn or Lp], and promote EHP collection, as depicted in Fig. 1(b). The 
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gradient interface between the P*/P and N/N* layers generate an extra electric field Ey and zb 
P N 


which reduces battery surface recombination and enhances EHP collection. In this 3D structure, 
EHPs are mainly generated in the ridges, and hence these areas contribute the most to the output 
power. Monte Carlo simulations and COMSOL Multiphysics were used to optimize the structural 
parameters, including the thickness of each doping region, doping concentration, ridge width, and 


ridge spacing, and predict the battery output performance. 


®3Ni source was selected as the beta source due to its long half-life (approximately 100 years), 
moderate decay energy (Eavg=17.4 keV, Emax=66.9 keV), and solid metal form, which allows easier 
and safer handling. SiC was selected as the converter semiconductor material because of its 
desirable properties including low leakage current density, higher radiation damage threshold, 
higher conversion efficiency, and excellent tolerance to harsh environments, i.e., extreme 
temperatures, wear, chemical exposure, and radiation [37-39]. Moreover, the development of SiC 
etching technology has enabled the fabrication of microgroove structures with high aspect ratio [39]. 
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Fig. 1. (a) Schematic 3D diagram showing part of the proposed battery and the distribution of the electron-hole pair 
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(EHP) generation rate in the innermost ridge with 0.8 um and d=1.2 um; (b) Structure principal diagram of SiC 


P*PNN* betavoltaic nuclear battery. Dimensions are not in scale. 

2.2. Methods 

2.2.1. Radiation-induced carrier generation in 3D diode structures 

A. Model for radiation-induced carrier generation rate in 3D diode structures 


The distribution of the EHP generation rate in a betavoltaic nuclear battery is governed by the 
energy deposition of beta particles in the converter, which significantly affects the output 
performance. In the conventional planar diode structure depicted in Fig. 2(b), the energy deposition 
[Eaep(x)] along the radiation transport depth [x] in bulk SiC was calculated via a Monte Carlo 
simulation with the Geant4 radiation transport toolkit, using a rectangular °Ni source with a full 
energy spectrum. The “Ni source is characterized by the specific activity of 5.68 Ci/g, 100% 


abundance, and density of 8.9 g/cm’, with isotropic emission of beta particles. SiC has the bandgap 


width of 3.26 eV, relative dielectric constant of 9.7, density of 3.21 g/cm’, and intrinsic carrier 
concentration of 7.4x10° cm®, derived using the widely employed formula [40]. G(x) is obtained 


and expressed as: 
G(x) = Aart = Go exp(—ax) (2) 


where A is the activity of the radioisotope, ¢ is the average energy needed to generate an EHP, 
commonly considered as 6.88 eV for SiC [41], and Go and a are the surface EHP generation rate 
and absorption coefficient, respectively. G(x) of Ni sources with varying thickness [f] is presented 
in Fig. 2(c), and the corresponding Go and a are derived by fitting G(x). The dependences of Go and 
a on the source thickness [f] are then fitted, as demonstrated in Fig. 2(d). Consequently, G(x) 
considering different source thicknesses is expressed as: G(x, t) = A,(1 — eit) @- (Are M2F+B, )x n 
the fitting parameters for which are provided in Table 1. Fig. 2(c) shows that the saturation thickness 
of the Ni source is 3 um, consistent with previous work [24, 42]. The penetration depth of B 
particles emitted from “Ni is 10.2 um in SiC, defined as the location where 99% of the total energy 


deposition occurs according to Alam [8]. 


To accurately model the energy deposition and EHP generation rate in a multi-groove 
betavoltaic nuclear battery, it is essential to consider the superposition contributions of all the 
isotope sources and ridges. This involves extending G(x,f) of the 2D structure to that of the 3D 
structure. Fig. 2 (a) characterizes the penetration distance distribution of B particles released by °Ni 
isotope sources, where D; denotes the jth layer converter (the jth ridge of the proposed battery) and 
S; indicates the jth “Ni layer. At a random location of the EHPs generated in converter D; (specified 
by reference site P), the distance between P and source Sj-1 is represented by [x]. Converter D; is 
exposed to [j-1] layers of isotope sources on its left, and the B particles emitted from the ith (i<j) 
source must traverse a total source thickness of [(j-i-1)¢] and a total converter thickness of [(j-i- 
1)d+x] to reach reference site P. Similarly, D; is exposed to [n+1-] layers of sources on its right, and 
the B particles emitted from the Ath (<= k <=n) source must penetrate a total source thickness of 
[(k-j)t] and a total converter thickness of [(A+1-/)d-x] to reach P. During this process, the energy of 
the beta particles decays exponentially with penetration distance, and the EHP generation rate of D; 
is given by: 


G(x, t) = Yi Yep (te OG Dats) vOG-- Oe} $ EKEN {Go(t)e tO (Oct 1 dv) e-v k-t) (3) 


where n is the layers of sources, Go(t) is the surface EHP generation rate, a(t) is the absorption 
coefficient of B-electron flux in SiC, and y(t) is the absorption coefficient of B-electron flux in Ni, 
which are acquired via the equation in Table 1 based on Monte Carlo code simulation. The EHP 
generation rate of the (n+1)th converter can be expressed as Gn+1(x, t) = G,(d — x,t), owing to 
the symmetry of the multi-groove structure. Table 1 shows that the R? (R-squared or coefficient of 


determination) values are greater than 0.99, indicating the excellent fitting performance of the model. 
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Table 1. Exponential fitting parameters and R-squared values for Go(t), a(t), and y(t). 


Equations Fitting values Ug 
Go(t) = Aet + B} Ay=-8.14« 108 cm°®-s!, u1=1.35 pm, B1=8.69x10!5 cm°- s”! 0.999 
a(t) = Azet + B, A= 3638.76 cm", u2=2.19 pm, B2=5216.14 cm! 0.992 
y(t) = Azet +B} A3= 8472.34 cm, u3=2.06 umt, B3= 15344.09 cm! 0.997 


15 
K d j-i-1)t+(j-i-l)d+x  (k-j)t+(k+1-j)d. (b) (eraa = 
t iHi BORA Jaa ti x Fit fuction G(x)=Gyexp(-ax) 
i 3 = =0.1 um 
i 1, e = 0.5 pm 
E 6 a f= 1.0 um 
D 2 v 3.0 um 
a & 4 <= 9.0 um 
4 1) 
0 
y 0 4 8 12 16 
© SiC mm Ni-63 e ß particle Penetration depth (um) 
15 3 Data from Geant4 
(d) 12 ai 2 aris x10'6 — G(x) by fitting formula (3) 
— Fit fuction Gp(x)=4 ,e""+B, b (e) t= 1.5m d=3um R? = 0.985 
10} = G, b | | oii.” 
=> 78 = 05+ 
bs 8 = ii S e È 72 18 24 30 36 42 
g 47 g Tta (f) t=1um d=3 um R? = 0.990 
2 6 — Fit function a(t) = A,e"?'+B, ~ 7 1.0F 
Ss = 4 3 R 
šl va 2 0. Lge ee eee ee 
4 oh 6 Š 0.0 L L L 1 == -i 
Oige 6 12 18 24 30 36 42 
ais 
2 iA (g) t= 1m d= 4m R? = 0.982 
i Me Gd i 3 7 1.0F 
0 2 4 6 8 10 È osl 
Source thickness,t (um) 5 00 , m , , 1 , , 
0 8 16 24 32 40 48 56 64 


Relative depth (um) 


Fig. 2. Schematics of (a) 3D and (b) 2D converters; (c) EHP generation rates along penetration depth and (d) surface 
EHP generation rate and absorption coefficient for different source thicknesses in 2D converters; EHP generation 
rates in 3D converter with (e) 1.5 um ridge spacing and 3 um ridge width, (f) 1 um ridge spacing and 3 um ridge 
width, (g) 1 um ridge spacing and 4 um ridge width. 


B. Validation of electron-hole pair (EHP) generation rate model in 3D diode structures 


Accurate prediction of the rate of EHP generation is essential for optimizing the design of 3D 
SiC-based betavoltaic nuclear batteries. Therefore, the validity and ability of the model to predict 
data from 3D multi-groove structures was demonstrated through comparison with original Geant4 
data [AE dep(x)/e]. Figs. 2 (e)-(g) exhibit the reliability of the proposed EHP generation rate model, 
expressed by formula (3), with highly consistent G(x,f) curves compared with the original Geant4 
data for both ridge spacings of 1.5 and 1 um, and ridge widths of 3 and 4 um. The high R? values of 
0.985, 0.990, and 0.982 for these curves indicate the accuracy of the model in calculating the EHP 


generation rates in 3D structures, rendering it a valuable tool for optimizing SiC-based betavoltaic 


20 


nuclear batteries for high performance. 


In the 2D converter, the EHP generation rate decreased exponentially with depth; whereas, the 
multi-groove 3D structure exhibits a unique distribution of higher EHP generation rates in the inner 
ridges and lower rates in the outermost ridges. The innermost ridges exhibit high EHP generation 
rates on the lateral surfaces, low rates in the middle, and a symmetrical distribution, matching the 
3D EHP distribution in the inner ridge shown in Fig. 1(a). These findings suggest that the proposed 
multi-groove 3D structure has the potential to significantly enhance the power output compared 


with the traditional 2D structure, particularly as the relative depth increases. 
2.2.2. Model for radiation-induced current in 3D diode structures 


Fig. 1(b) illustrates that radiation-induced EHPs generated within the depletion region can be 
collected with 100% efficiency, whereas those generated outside the depletion region can only be 
collected after diffusion to the PN junction boundary, the P*/P interface, or the N/N* interface. The 
CE(y) was calculated using the equation [43]: 


ay) 


CE(y) = 1 — tanh 7 


(4) 


where d(y) represents the distance from the depletion region boundary or the interfaces and is set to 
zero inside the depletion region. The electron and hole minority carrier diffusion lengths L are 
expressed as Ini, Ln, Lp, and Lp, in regions Pt, P, N, and N+, respectively. The radiation-induced 


current density Jp can be expressed as: 


Jr = (LEM g J. Gj (x,t) dx fy" CEC) dy Jf.” dzy/S (5) 
where q is the electron charge; H is the ridge width [d]; and H2 and S are the rectangular side length 
(1 cm) and area of the device (1 cm’), respectively. 

According to the drift-diffusion theory, nonequilibrium carriers generated within the depletion 
region and the neutral region outside the depletion region boundary within a minority diffusion 
length can be collected, thus contributing to the current density (Jr). The effective charge collection 
region (ECR) length, represented by H_ECR= (Wat Lni+2Ln+2Lp+ Lp1), determines Jr and can be 
maximized by increasing Wa, Ini, Ln, Lp, and Lp1. Lower doping concentrations increase Wa and L, 
leading to higher EHP collection, as demonstrated in Figs. 3(a) and (b), where Wa decreases from 
7.58 um to 30 nm and Ly decreases rapidly from 77.34 to 11.06 um with increasing doping 
concentration from 10! to 10!° cm”. Lp decreases more gradually over the same doping range. The 
minimum doping concentration in the P- and N-regions was 10'4 cm? due to our facility’s capacity 
to process low doping in SiC materials; whereas, the minimum doping concentration of the heavily 
doped P*- and N*-regions required Na =10!° cm? and Np =10!8 cm”, respectively, to reduce the 
ohmic contact. To maximize H_ECR and Jp, the maximum values of Wa, Ini, Ln, Lp, and Lp: should 


be adopted, as listed in Table 2 and detailed in Supplementary Materials S1 and S2. 
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Fig. 3. (a) Depletion region width and (b) minority diffusion length with doping concentration. 


Table 2. Maximum depletion region and minority carrier diffusion length in each region of the converter. 


Wa_max (um) Lni_max (um) Ln_max (um) Lp_max (um) Lp1_max (um) 


7.58 11.06 77.34 11.44 6.15 
2.2.3. Battery output characteristic model and simulation 


The electron (hole) concentration [n (p)] inside the semiconductor converter device is 


governed by the carrier continuity equation as follows: 


ên at ee 6 
ot q Jn n (6) 
ae +G—R (7) 
ot q p p 


where jn (jp) is the electron (hole) current density, G is the EHP generation rate derived from Eq (3), 
and R, (Rp) is the electron (hole) recombination rate. jn (jp) can be described by the drift and diffusion 


processes [44] as: 
Jn = NAHE + qD,Vn (8) 
Jp = PAHpE — qDpVp (9) 
where Æ is the sum of the external and internal electric fields generated by the diffusion of EHPs. 
These relationships are governed by the Poisson equation: 
gf a a =e) (10) 
Ox E0€Er 
where p is the charge density, €o is the vacuum permittivity, €, is the relative permittivity of the 
semiconductor, and Np (Na) is the concentration of donor (acceptor). 
A. Current-voltage characteristics of ideal PN junction diodes 
The current-voltage (J-V) characteristics of an ideal long PN junction diode can be derived 
assuming that the external voltage drops entirely in the depletion region, the diode operates under a 
low injection level (i.e., the excess carrier concentration is much smaller than the equilibrium 
majority carrier concentration), there is no recombination or generation current in the space-charge 
region, the semiconductor is nondegenerate, and the length of the P- (N-) region is much larger than 


the diffusion length[45, 46]. 


Under these assumptions, the J-V relationship of the neutral areas (E=0) of the P- and N- 
regions can be obtained by solving the continuity equations (6)—(7) and current equations (8)-(9). 
Ignoring the recombination and generation currents in the space-charge region, the J-V 


characteristics of an ideal diode are given by: 


J =Jsc — Up + Jn) = Isc — Jo [exp (=) -1 (11) 


where Jsc is the short current, equals to radiation-induced current Jp; k is Boltzmann’s constant; T 
is the absolute temperature (300 K in this work); V is the bias voltage; and Jo is the leakage current 


density of the PN junction. 


The open-circuit voltage is derived as J=0 A. 
kT (Jsc (12) 
Voc = — lIn (= + 1) 
c qg Jo 
Then, the maximum output power Pm can be described as: 

Pn = JscVocFF (13) 
where FF is the fill factor, which can be derived using the open-circuit voltage [13]. Jo is determined 
by the diffusion efficiency and diffusion length of the minority carriers, and the doping 
concentrations, shown in S3 of the Supplementary Materials. Using this ideal diode model, the 
battery performance can be predicted quickly, and the optimal battery structure can be determined 


based on the J-V characteristic numerical model. 
B. Current-voltage characteristics calculation using COMSOL Multiphysics 


The J-V characteristics of real diodes differ from those of ideal diodes owing to the various 
assumptions made for ideal PN junction diodes. Therefore, an accurate method is essential to 
simulate the J-V characteristics of real diodes, considering the generation, recombination, and drift 


of charge carriers, as well as the real characteristics of the semiconductor material. 


COMSOL Multiphysics is a powerful tool for simulating the current-voltage characteristics of 
realistic situations by solving partial differential equations that incorporate real physical phenomena. 
The simulation utilizes various physical models, including the Monte Carlo simulation for 
calculating the EHP generation rate in 3D diode structures (Equation (3)), the Shockley-Read-Hall 
model for trap-assisted recombination (Equation (14)) [45], and the low-field mobility model for 
determining the minority carrier mobility (Equations $4-S5 in the Supplementary Materials). 
Additionally, the carrier lifetime model is described in Equations S6—S7 in the Supplementary 


Materials. 


CyCaN,(np — n?) 


np — n? 


R= 


(14) 


Here, Cy (Cp) is the electron (hole) capture coefficient, n; is the intrinsic carrier concentration, Tn (tp) 
is the electron (hole) minority lifetime given by Equations S6—S7 of the Supplementary Material, 


E; is the recombination center (defect) level, and Æi is the intrinsic Fermi level. The defect energy 


: 8 1 1 
level is set as E, = E;, and the defect density is Ny = z ar 
nTn ptp 


Solving partial differential equations (6)-(10), we obtained important information on the 
electric field, carrier recombination, and electron (hole) current density, enabling prediction of 
crucial electrical parameters, including the short-circuit current [Jsc], open-circuit voltage [Voc], 
and output power [Pm]. Ultimately, with this approach, the behavior of diodes can be accurately 
modeled, which is vital for optimizing their performance and ensuring that they meet the 


requirements of various real-world applications. 
3. Results and discussion 
3.1. Performance calculated using ideal diode model 


To optimize the structure of the proposed battery and maximize its output power density, a 
parametric sweep was conducted in the numerical model to adjust variables including the single- 
source thickness (ridge spacing, f); single-converter thickness (ridge width, d); thicknesses of the 
P*-, P-, N-, and N*-regions (H_P*, H_P, H_N, and H_N*); acceptor concentration of P-region (Na); 
and donor concentration of N-region (Na). Possible values of ¢ and d are in the range of 0.1-10 um 
with a 0.1 um step size for each. It is worth noting that the source thickness is fixed to the saturation 
thickness of 3 um as ¢ is greater than 6 um, and evenly distributed on the inner surface of the 
microgroove. The feasible ranges for H_P and H N are 1 to 250 um and 1—60 um, respectively, 
while Na and Na can range from 1x 10!* to 7.94 10!8 cm? and 1x10! to 7.94x10! cm >, respectively. 
The values for H_P* and H_N* correspond to 10 and 6 um, which are close to the minority carrier 
diffusion length in the P*-region and N*-region. Additionally, the heavily doped P*-region and N*- 
region are assigned acceptor concentration and donor concentration values of 1 x10!° and 1 x10!8 
cm, respectively. 

3.1.1 Optimizing ridge spacing and width 

The output power of a betavoltaic nuclear battery depends on the amount of beta particle 
energy deposited in the converters and the efficient collection of radiation-induced EHPs. The 
coupling between the radiation source and the device is crucial for enhancing the output power. To 


illustrate the impacts of ¢ and d on the output performance of the betavoltaic battery, the 3D surface 


contours of the short-circuit current density (/sc), open-circuit voltage (Voc), and maximum output 


power density (Pm) were plotted, as shown in Fig. 4. 


As depicted in the 3D surface wireframe perpendicular to the d direction in Fig. 4(a), Jsc 
initially increases rapidly, but subsequently decreases as ¢ increases, reaching a peak value with ¢ of 
0.1—2.2 um. Additionally, the dependence of Jsc on d exhibits a similar trend, reaching its peak with 
d of 0.2-4.0 um. The bottom-projected contour shows that the optimal values of Jsc are achieved 
when d is in the range of 0.2-3 um and t ranges from 0.1 to 2 um. The maximum Jsc of 8.57 A/cm? 


is achieved by combining 0.8 um and d=1.2 um. 


Fig. 4(b) demonstrates that Voc initially increases as ¢ increases and reaches saturation at =2 
um, while it decreases with an increase in d. Nevertheless, overall, Voc is not highly sensitive to 
variations in ¢ and d, with a range of 2.36—2.46 V. Due to the minor variation in Voc, Pm is primarily 
determined by Jsc. The relationship between Pm and variations in ¢ and d exhibit similar trend to 
that of Jsc, with a rapid initial increase and a subsequent slow decrease with increasing ¢ or d, as 
illustrated in Fig. 4(c). The maximum Pm value of 19.73 uW/cm? is attained at 0.8 um and d=1.2 


um. 


The factors affecting the output power mainly include A, ys, (1-r), Q, Voc, and FF. With 
varying t and d, A, 75, (l-r), Voc, and FF exhibit relative increases of 9844.28%, 578.19%, 1956.94%, 
4.14% and 0.21%, respectively, while O remains unchanged at 45.64%. Therefore, the key factor 
affecting the output power is the coupling of A, ys, and r. For the detailed calculation process, please 


refer to Supplementary Material S4. 


The input power [Pin] can be calculated by combining A, ys, and r, as Pin= AF avens(1-r). Figs. 
4(d)-(e) present the maximum output power [Pm] and its corresponding optimized ridge width [d] 
of 0.2—4.0 um, for different source thicknesses [t]. We found that Pm is determined by the coupling 
of A, ys and r, i.e., Pin. The optimized value of d gradually saturates with increasing t, consistent 
with the phenomenon of energy deposition saturation in the converter with increasing d. 
Additionally, the overall conversion efficiency [7 tot] increases and then decreases with t, reaching a 
saturated value of approximately 1.5%, and mot corresponding to maximum Pin and Pm is 4.58%. 

In the proposed battery, smaller ridge spacing and ridge width can improve the source activity 
and reduce the self-absorption effect. However, excessively thin converters may result in a lower 
(1-r) and do not match the particle penetration depth, leading to a reduced Pin. Therefore, a tradeoff 


between ¢ and d is necessary to maximize the power density. 
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Fig. 4. Effects of the ridge spacing [7] and ridge width [d] on the (a) Jsc, (b) Voc, (c) Pm; the dependences of (d) Pm, 
Pin, (€) d, and Ntot on t. Na and Na are both set to 1 x10'4cm™; H_Pand H_N are set to 158 um and 26 um, respectively, 
approximately equal to (2Ln_ max + Wa_max /2) and (2Lp_max + Wa_max /2). 


3.1.2 Optimizing widths of P-region and N-region 

After optimizing ¢ and d, we investigated the dependence of Jsc, Voc, and Pm on the widths 
of the P- and N-regions for betavoltaic nuclear batteries, as shown in Figs. 4(a)-(c). Increasing H_P 
results in an initial increase and subsequent saturation of Jsc, Voc, and Pm. This trend can be 
attributed to the increase in the radioisotope source activity with H_P, which generates more EHPs 
in the ECR. However, when H_P exceeds the ECR length, the performance metrics reach saturation. 
Specifically, Jsc, Voc, and Pm saturate when H P reaches 160 um, and an additional increase of 10 
um results in a negligible increase of less than 1%. Similarly, Jsc, Voc, and Pm present a similar 
trend with H_N, showing a gradual rise and saturation at H_N =24 um. This phenomenon can be 
attributed to the longer minority carrier diffusion length (Zn=77.34 um) of the P-region compared 
with that of the N-region (Zp=11.44 um), and their saturation values of H_P and H N are around 
twice Ln and Lp, respectively. Therefore, H_source should not exceed 200 um, where H_P* and 
H _N* are 6 and 10 um, respectively, and H_P and H_N should not exceed their saturation values, 
corresponding to 160 and 24 um. 

Fig. 5(d) shows that A linearly increases with H_P from 0.51 to 6.59 Ci (by 1204.00%), while 


Q first increases and then decreases with the relative change rate of 164.92%. The optimization of 


H P reveals that factor A is the dominant factor affecting the output power, and factor Q is secondary. 
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Although A and Q show different trends with changes in H_P, their products, A-O and Pm, exhibit 
consistent trends, further confirming that A and Q are the key factors affecting Pm, as shown in Fig. 
5(e). However, Voc and FF increase minimally by 1.81% and 0.09%, while 7, and (l-r) remain 
unchanged at 61.53% and 48.77%, respectively. Please refer to Supplementary Material S4 for 


detailed calculations. 


As depicted in Fig. 5(f), for the structure with Æ source = 200 um, Jsc and Pm initially 
increase with H_P and then decline instead of continuously increasing. This behavior is due to the 
contribution of H_N to Jsc, and Pm outweighs that of H_P when H P approaches its saturation 
thickness. Therefore, when optimizing the widths of the P- and N-regions, the balance between H_P 
and H_N should be considered, as the combination of H_P =156 um and H_N =28 um achieves the 
maximum Pm of 19.74 u.W/cm*. A longer P-region is more suitable for a betavoltaic nuclear battery 
with a P'PNN* junction structure because it is more conducive to enhancing the collection of EHPs 
owing to the larger minority carrier lifetime and mobility compared with the N-region, as shown in 


Fig. S2 in the Supplementary Material. 
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Fig. 5. Effects of H_P and H_N on (a) Jsc, (b) Voc, and (c) Pm with varying H_source; dependences of (d) A, Q, 
(e) Pm, and A-O on H P with H N fixed at 24 um; (f) dependence of Pm on H_P with H_ source = 200 um. Na and 
Na are both set to 1 x10'4 cm, t=0.8 um, and d=1.2 um. 


3.1.3 Optimizing doping concentration in P- and N-regions 
Based on the analysis in section 2.2.2, the doping concentration plays a crucial role in 


determining the depletion width and minority carrier diffusion length of the betavoltaic nuclear 


Voc (V) 


battery, which significantly affects the collection efficiency of the EHPs and the output performance 
of the battery. 


Figs. 6(a)—5(c) illustrate the effects of Na and Nz on Jsc, Voc, and Pm for the proposed battery. 
As shown in Fig. 6(a), Jsc increases with a decrease in Nz owing to its beneficial effect of expanding 
the minority carrier diffusion length and depletion region width to promote EHP collection, as 
depicted in Fig. 3. The variation in Jsc with Na was small compared with N, because Lp is much 


smaller than Ln, resulting in less changes in the collection efficiency of EHPs. 


Fig. 6(b) shows that Voc initially increased rapidly but subsequently decreased slowly as Na 
increased, as depicted in the 3D surface wireframe perpendicular to the Na direction. The 
dependence of Voc on Na exhibits a similar trend. The maximum Voc of 2.66 V is obtained by 
combining N, =3.1610!8 cm? and Na=7.94x10!7 cm”. The bottom-projected contour indicates that 
the optimal values of Voc are achieved at higher doping concentrations, which is attributed to the 
reduction in the leakage current [Jo] owing to the higher doping concentration, as demonstrated in 


Fig. S3 in Supplementary Material. 


Although Jsc and Voc exhibit opposite dependencies on doping concentration, Pm varies in 
the same way as Jsc because the variation range of Voc with doping concentration is small, of 2.44— 
2.66 V with a relative increase of 8.98%. The optimal doping concentration combination is Na 


=1x10'4 cm? and Na = 1x10!4 cm”, yielding the maximum output power density of 19.74 uW/cm?. 


As shown in Fig. 6(d), low doping increased Q from 13.48% to 45.65% owing to the widened 
depletion region and diffusion length, leading to longer ECR length [H_ECR] and improved EHP 
collection, and resulting in a higher power density. H_ ECR and Pm follow a similar trend as Q, 
confirming that low doping enhances the power density by increasing H_ECR for EHP collection. 
However, at the highest doping concentration, H_ECR and Pm are minimum, of 53 um and 6.37 
uW/cm?, respectively, which are only slightly better than the power of 5.80 uW/cm? with the 
H source of 53 um. Pm only marginally improves even with a nearly four-fold increase in the source 
activity [A], indicating the critical role of H_ECR in the design of the 3D battery and that H_source 
should not exceed H_ECR. 


Additionally, FF ranges from 94.18% to 94.57%, with a tiny relative increase of 0.42%, while 
A, ys, and (l-r) remain constant at 4.04 Ci, 61.53%, and 48.77%, respectively, resulting in constant 
Pin of 129.16 uW/cm?. Therefore, Q is a significant factor affecting the output power density, and a 
low doping concentration in the P- and N-regions is recommended to enhance it, hence maximizing 
the short current and output power density. Detailed calculations are provided in Supplementary 


Material S4. 


202310.00022v1 


chinaXiv 


(a) 


Isc (uA/cm?) 
$ 
x 
`” 
a =¥— Ng=1e14 /cm? 
2 — 200 9 g=1e14 /cm 
Pm (uW/cm~) E 460 -4- Ng=1e16 /cm3 
= Š 120 
E =! gob te Ng=1e14 Jem? 
š =å- Ng=1e16 /cm? 
a& 
g —+ Na=1e14 /cm? 
E -@ N4a=1e16 /cm? 
Š 
= 
3 
E 
a 
5 
Bio SY 14 15 16 17 18 19 


4a) 
“Mo (cm 53) 


Fig. 6. Influence of Na and Na on (a) Jsc, (b) Voc, and (c) Pm; effect of Na on (d) H_ECR, Q, and (e) Pm. H_P = 156 
um, H N=28 um, 0.8 um, and d=1.2 um. 


3.1.4 Critical parameters of device structure 


To achieve the best device performance, we comprehensively optimized the geometric 
dimensions (optimization procedure #1), doping concentration (optimization procedure #2), and 
width of each doping layer (optimization procedure #3) of the semiconductor materials. The final 
optimized parameters for the battery are Na = Na =1 x10!4 cm, H_P=156 um, H_N=28 um, 0.8 
um, and d=1.2 um. 

In Optimization procedure #1, A, ys, and (l-r) exhibited significant relative increases, while 
the relative increases in Q, Voc, and FF were negligible, as shown in Fig. 7(a). Therefore, the key 
factor affecting the output power is the coupling of A, ys, and (1-r), and a tradeoff between them is 
required to maximize the power density, as depicted by their optimal parameters corresponding to 
the maximum output power. 

In Optimization procedure #2, A was identified as the dominant factor affecting the output 
power, with the significant relative increase of 1204.00%, and Q was a secondary factor with the 
relative increase of 164.92%, as depicted in Fig. 7(b). Increasing H_P and H_N can enhance the 
source activity [A] and ECR length. However, increasing H_P and H_N beyond the ECR leads to 
Q, as EHPs outside the ECR become difficult to collect. The tradeoff between A increasing and Q 


dropping results in saturation of the output power. 


In Optimization procedure #3, we found that Q is the key factor affecting the output power; 


whereas, the other factors have negligible relative increases, as shown in Fig. 7(c). To optimize the 
output power density, low doping concentrations should be used in both P- and N-regions. This 
enables larger diffusion lengths and wider ECRs, resulting in a higher collection efficiency of EHPs, 
particularly through an extended H_ECR. 


Among these parameters, the coupling of A, ys, and (1-r) has the greatest impact on the output 
performance, which is related to the activity, self-absorption, and backscattering of the radioactive 
source. The second most influential parameter is Q, which depends on the depletion region width 
and diffusion length controlled by the doping concentration, as well as the ECR length governed by 
H_P, H_N, and the doping concentration. To maximize the output power density, use of thinner 
sources and converters, lower doping concentrations, and larger H_P, of approximately twice the 


diffusion length, are recommended. 


To validate the numerical model established in this study, we performed calculations on 
planar batteries with the same geometric dimensions and semiconductor parameters as those in 
references [24] and [16]; good agreement was apparent regarding the output power density, as 
shown in Table 3. It is worth noting that the power density of battery #1 was converted based on the 
100% abundance of Ni source, in contrast to the 20% abundance mentioned in reference [24]. The 
difference in the results between reference [16] and this study (#5 battery) is relatively large because 
the previous work did not consider the energy loss caused by the collection efficiency Q, leading to 
an overestimation of Pm. Using our method, we can further provide suggestions for optimizing the 
batteries in references [24] and [16] by adjusting Na, Na, H_P and H_N. For reference [24], setting 
H P, H N, Na, and Na to 4.9 um, 19.9 um, 1x10! cm”, and 3.98x10'4 cm”, respectively, the 
maximum output power can be increased to 296.5 nW/cm”, representing a 12.6% improvement as 
depicted by battery #3. For reference [16], by combining H_P =0.1 um, Na and Na at 4.7x10'3 em? 
and 3x10!6 cm”, respectively, the output power can be increased to 378 nW/cm?, with a 10.5% 


improvement as depicted by battery #6. 


Finally, we conducted a comparative analysis of the P*PNN* and PN structures, assessing 
their respective output power densities in the 2D (#7 and #8) and 3D (#9 and #10) batteries. To 
achieve this, we substituted the P* and N* layers in the P*PNN* structure with P and N layers of the 
same size to obtain the PN structure, as shown in Table 3. Our results indicate that the P*PNN* 
structure outperforms the PN structure, with a 10% increase in Pm in the 2D configuration and a 39% 
increase in the 3D configuration. These findings highlight the superior performance of the P*PNN* 


structure, particularly in 3D configurations. 


In this study, the diffusion lengths were calculated through commonly used equations [17, 
42], and the minority carrier lifetimes (tn and tp) adopted in these calculations (shown in 
Supplementary Material S4) fall within the range of experimental values (tn of 0.9—10 us and tp of 
0.05—2.1 us) [47-50]. Moreover, experimental diffusion lengths for SiC were reported to be 30—100 


202310.00022v1 


chinaXiv 


um [49, 51], which provides further validation for the calculated electron and hole diffusion lengths 


presented in this article. 


(Table 3. Comparison between our battery model results and previous works; optimization suggestions for the 


batteries in references; comparative analysis between the P*PNN* structure and the PN structure. 


Battery Source thickness Converter Doping 
Type Pm (nW/cm?) Ref 
ID (um) thickness (um) concentration(cm™) 
#1 2D-P*PNN* 1 269.1 [24] 
P*\P\N\N™: P*\P\N\N*: 
m 2D-P*PNN* 1 0.1\23.9\0.9\0.1 —-10!9\10!\10!9\1019 263.3 Qur work 
for validating 
P*\P\N\N™: P*\P\N\N*: Our work 
#3 2D-P*PNN* 1 296.5 
0.1\4.9\19.9\0.1  10!9\10'4\3.98x10!4\10!° for optimizing 
#4 PN 2 360 [16] 
P\N: P\N: 
13 16 Our work 
#5 PN 2 8\15 4.7x10'*\3x10 342 
for validating 
Our work 
#6 PN 2 0.1\22.9 4.7x10!\3x10!6 378 
for optimizing 
P*\P\N\N™: P*\P\N\N*: 
#7 2D-P*PNN* 2 394 
0.1\4.9\19.9\0.1  10!9\10'4\3.98x10!4\10!° 
#8 2D- PN 2 P\N: 5\20 P\N: 10!4\3.98x10!4 358 Our work 
P*\P\N\N': P*\P\N\N™: for comparing 
#9 3D-P*PNN* 0.8 19.74 W/cm? 
10\156\28\6 10!\10!*\10'4\10!8 
#10 3D-PN 0.8 P\N: 166\34 P\N:10!4\10"4 14.22 wW/cm? 
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Fig. 7. Optimal parameters and their relative increase considering different optimization procedures: (a) Optimizing 
t and d; (b) Optimizing H_P and H_N; (c) Optimizing Na and Na. Relative increase is calculated as the difference 
between the maximum and minimum values divided by the minimum value. The bars show the minimum values 


only, indicating that the maximum value is equal to the minimum. 
3.2. Performance simulation using COMSOL Multiphysics 


The largest discrepancies between the ideal and practical performances of betavoltaic nuclear 
battery are the collection efficiency Q and Voc FF/e [10]. To better understand the discrepancies 
between the ideal and practical performance of betavoltaic nuclear batteries, it is necessary to 
conduct a detailed analysis of the specific differences caused by Q, Voc, and FF, as well as their 
underlying reasons. COMSOL Multiphysics was employed to simulate the practical performance of 
betavoltaic batteries. Most SiC material properties were imported from the COMSOL library; some 
properties that were not available in the library, such as minority carrier mobility and minority 


carrier lifetime, were manually added based on the literature [17]. 


In the simulation process, a single converter was modeled with dimensions of 1 um x 1.2 um 
x 200 um. A user-controlled mesh of four different sizes (#a, #b, #c, and #d) was defined to improve 
both the accuracy and computation time, with maximum element sizes of 50, 100, 200, and 500 nm. 
The minimum element sizes were set to 1/10 of the corresponding maximum element sizes. 
Although the use of finer mesh sizes resulted in slightly larger Pm values, the differences are 
negligible. Specifically, Pm values obtained from mesh #a are only 0.064%, 0.121%, and 0.149% 
higher than those obtained from meshes #b, #c, and #d, respectively. Based on these results, the #d 
mesh was chosen for computation, with boundary elements set at the maximum element size of 0.1 
and minimum element size of 0.02. The maximum element growth rate was set to 1.1 with the 


curvature factor of 0.25, and the resolution of the narrow regions was specified as 1. 


Fig. 8(a) displays the relationship between Pm and H_P obtained through numerical modeling 
and COMSOL simulations. The variations in Pm with H_P from these two methods are in excellent 
agreement, reaching their maximum values at H_P of 156 and 164 um, respectively, with values of 
19.74 and 18.69 uW/cm?, differing by only 5.62%. However, there were some differences in Pm at 
lower values of H_P, which increased as H_P decreased. These differences arise from Q, Voc, and 
FF, as shown in Figs. 8(a) and (b), respectively. The numerical model relies on empirical formulas 
(4) and (12) for O and Voc, resulting in a lower Q and higher Voc than the COMSOL simulation, 
which calculates Q by subtracting the Shockley-Reed-Hall recombination rate from the EHP 
generation rate (shown in Fig. 8(d)) and extracting Voc by finding the voltage value when the current 
is nearly zero. In addition, the numerical model produces a higher FF with little variation, while the 
COMSOL simulation calculates FF by dividing Pm by Voc and Jsc, which are sensitive to the 


distribution of Voc. 


In addition, Pm initially increases and then decreases with increasing H_P, reaching a 


maximum value at H P of approximately 2Zn. This can be explained by the suppression of carrier 


recombination owing to the electric field distribution. The electric field distribution intensity and 
range at H_P=164 um are much greater than those at H_P=4 um, resulting in a considerably lower 
carrier recombination rate, leading to a larger Q and hence a larger Pm, shown in Fig. 8(c). This 
indicates that the internal mechanism affecting the power increase is the electric field distribution, 


which governs the carrier transport and collection characteristics. 


The results demonstrate that the proposed battery has a significantly improved output 
performance compared with conventional planar batteries reported in previous studies [13-16, 19, 
23, 24], with a maximum output power density approximately 50 times higher. If multiple batteries 
are stacked to a total thickness of around 10 mm, the proposed battery would provide an output 
power of approximately 1 mW, which can satisfy the power requirements of MEMS with 
dimensions less than 10 mm and power consumption of 1-100 uW. We provide further 
recommendations for designing 3D-groove betavoltaic nuclear batteries regarding 3D etching 


techniques on SiC, as shown in Table S2 of the Supplementary Material. 
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Fig. 8. Discrepancies between the ideal and practical performances simulated by COMSOL: (a) Dependences of Pm, 
Q, (b) Voc, and FF on H_P; (c) Distribution of electric field and carrier recombination rate; (d) Dependences of (G- 
R) distribution on H_P and d. (G - R) represents the difference between the generation rate and the recombination 


rate of EHPs. 


4. Conclusions 


In summary, this paper presents a novel °Ni-SiC-based P*PNN* 3D structure with a multi-groove 
design that eliminates the need for preparing a PN junction on the inner surface of the microstructure 
and improves the performance of betavoltaic nuclear batteries. The fully coupled model developed 
in this study considers various factors, such as B-particle generation, self-absorption, backscattering, 
energy deposition, as well as radiation-induced carrier generation and drift-diffusion, and thus 
provides a valuable tool for efficient design and development of betavoltaic nuclear batteries with 
complex 3D structures. The epitaxially grown graded P/N layer significantly enhances H_ECR, 
promoting radioisotope source activity and carrier collection efficiency, producing the maximum 
output power density of 19.74 W/cm’, with relatively thin radioisotope sources and converters 
(t=0.8 um and d=1.2 um), lightly doped P- and N- regions (Na=Na=10'4 cm), and longer P-region 
widths (H_P=156 um). The analysis of carrier transport and collection characteristics using the 


COMSOL Multiphysics code provides insights into the internal mechanism of the power increase 


and clarifies the discrepancies between the ideal and simulated performances of betavoltaic nuclear 
batteries. However, the diffusion length is susceptible to process variations and a short diffusion 
length may reduce the advantages of the proposed P'PNN* 3D structure. In conclusion, the proposed 
3D structure with a multi-groove design combined with a fully coupled model and optimization 
methods presents a promising approach for designing and optimizing high-performance betavoltaic 
nuclear batteries. It is worth noting that the importance of H_ECR cannot be overstated by 
increasing the output power, as it directly affects both radioisotope source loading and charge 
collection efficiency. Furthermore, the 3D structure proposed herein is expected to be well-suited 


for narrow-bandgap semiconductor materials with ultralong diffusion lengths. 


See the supplementary material for details of our calculation procedure and recommendations on 
the device design of betavoltaic nuclear batteries considering three-dimensional etching techniques 


on SiC. 
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